New Dirac points and multiple Landau level crossings in biased trilayer graphene 
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Recently a new high-mobility Dirac material, trilayer graphene, was realized experimentally. The 
band structure of AB A-stacked trilayer graphene consists of a monolayer-like and a bilayer-like 
pairs of bands. Here we study electronic properties of ABA-sta,cked trilayer graphene biased by a 
perpendicular electric field. We find that the combination of the bias and trigonal warping gives 
rise to a set of new Dirac points: in each valley, seven species of Dirac fermions with small masses of 
order of a few meV emerge. The positions and masses of the emergent Dirac fermions are tunable 
by bias, and one group of Dirac fermions becomes massless at a certain bias value. Therefore, 
in contrast to bilayer graphene, the conductivity at the neutrality point is expected to show non- 
monotonic behavior, becoming of the order of a few e^ /h when some Dirac masses vanish. Further, 
we analyze the evolution of Landau level spectrum as a function of bias. Emergence of new Dirac 
points in the band structure translates into new three-fold-degenerate groups of Landau levels. This 
leads to an anomalous quantum Hall effect, in which some quantum Hall steps have a height of 
3e^ /h. At an intermediate bias, the degeneracies of all Landau levels get lifted, and in this regime 
all quantum Hall plateaus are spaced by /h. Finally, we show that the pattern of Landau level 
crossings is very sensitive to certain band structure parameters, and can therefore provide a useful 
tool for determining their precise values. 



PACS numbers: 73.22.Pr, 72.80.Vp, 71.70.Di 
I. INTRODUCTION 

Single- and multi-layer graphene materials are gapless 
2D semimetals with unusual band structure, in which 
low-energy excitations have non-trivial Berry's phases.^ 
In particular, low-energy excitations in graphene are two 
species (valleys) of massless Dirac-like fermions with lin- 
ear dispersion and Berry's phase tt. The excitations in 
bilayer graphene are chiral fermions with aparabolic dis- 
persion relation and Berry's phase of 2tt."' The unusual 
nature of excitations underlies fundamental phenomena 
found in these materials, including Klein tunneling, weak 
anti-localization, anomalous quantum Hall effects, as well 
as novel symmetry-broken states in high magnetic field.^H 

One important distinct feature of bilayer graphene (as 
compared to monolayer graphene and other 2D electron 
systems) is a unique tunability of its band structure. By 
applying a perpendicular electric fielcP^ which induces 
an asymmetry between the top and bottom layer, it is 
possible to induce a band gap of up to 250 meV.I^ Thus, 
bilayer graphene can be turned into a semiconductor with 
a widely tunable band gap, which makes it attractive for 
various device applications, and provides a way to study 
chiral carriers in new regimes. 

Very recently, a new carbon-based semimetaL trilayer 
graphene, has been realized experimentally.'^"^ Experi- 
mentally studied trilayer graphene samples had a high 
mobility, and exhibited quantum Hall effect (QHE). In- 
spired by these experiments, in this paper we study the 
effect of perpendicular electric field on the electronic 
properties of trilayer graphene. First, we aim to under- 



stand how the band structure of trilayer graphene trans- 
forms under bias, and to explore new possibilities offered 
by its tunability. Second, we would like to develop an un- 
derstanding of Landau level (LL) evolution and LL cross- 
ings under bias. Given that band structure of trilayer 
graphene is quite complicated, it is necessary to first un- 
derstand single-particle effects in LL evolution, before ad- 
dressing interaction phenomena, including quantum Hall 
ferromagnetism and fractional QHE. 

Depending on the stacking {ABA or ABC), the band 
structure of trilayer graphene can be very differentPsE^l^ 
Throughout the paper, we will focus on t he case of 
Bernal-stacked (ABA) trilayer graphene.!^ In the ab- 
sence of bias, band structure of Ai?A-trilayer graphene 
can be viewed as a combination of independent overlap- 
ping monolayer-like (Be rry's phase tt) and bilayer-like 
bands (Berry's phase 27r) piMHI Both bands are gapped, 
in contrast to gapless monolayer and bilayer graphene, 
however, the monolayer-like and bilayer-like gaps are 
overlapping such that there is no band gap. In contrast, 
ABC-stacked trilayer graphene has a Dira c poin ts with 
cubic band touching and Berry's phase StTj^^HII gpht by 
different symmetry-breaking terms. 

Previously, it was noted that perpendicular electric 
field completely transforms the band structure of trilayer 
graphene, hybridizing the monolayer-like and bilayer-like 
bands, ^ and generally opening only a small band gap. 
Here we explore the transformation of the band struc- 
ture in detail. Most interestingly, we find that as a result 
of trigonal warping, at bias ~ 0.5 V/nm,^^ seven new 
Dirac points (DP) emerge in each valley. One of the DPs 
is situated at the (K^) point in the Brillouin zone; 
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the remaining six consist of two groups of three, related 
to each other by 27r/3 rotations. The emergent Dirac 
fermions are massive, with both their mass and velocity 
being tunable by an electric field. At electric field ^ 0.9 
V/nm,l^ the band gap closes for three DPs. For higher 
bias the mass of these DPs changes sign, and overall gap 
increases monotonically. 

We analyze the evolution of LLs under bias, finding a 
very complex pattern of levels crossings and splittings. 
We illustrate the extreme sensitivity of this pattern with 
respect to band structure parameters. Even small varia- 
tions of certain parameters change the pattern of cross- 
ings at small bias. Therefore, the Landau level crossings 
should provide a useful tool for precise determination of 
the band structure parameters, which are still a subject 
of debate.^ The transformation of the band structure 
at high bias leads to an emergent three-fold degeneracy of 
some Landau levels in magnetic field. It is worth noting 
that the evolution of LLs with bias has been previously 
discussed by McCann and Koshino in Ref. [131 who, how- 
ever, considered an idealized model that neglected several 
smaller hopping parameters (see also Ref. jlQ). Here we 
consider the full tight-binding model,!^ finding that the 
smaller parameters play an important role in determin- 
ing the pattern of LL crossings. 

Our results also imply that single-particle physics in 
ABA trilayer graphene is much more complicated as com- 
pared to other chiral materials, such as monolayer and 
bilayer graphene, and it should be understood before ad- 
dressing the interaction effects. On the positive side, the 
tunability of the band structure may prove advantageous 
in studying and tuning correlated states, including frac- 
tional quantum Hall states , as discussed in the end 
of the paper. 

The rest of the paper is organized as follows: in the 
next Section we introduce the tight-binding description 
for the AB A-staxiied trilayer graphene and review de- 
coupling of its band structure into monolayer-like and 
bilayer-like block^. In Section HI we discuss the evolu- 
tion of the band structure as a function of perpendicular 
electric field. Landau level spectrum as a function of 
magnetic field and electric field is studied in Section [rv] 
Finally, in Section [Vj we provide a brief summary of our 
results and present a broader outlook. 



II. MODEL AND REVIEW OF UNBIASED 
CASE 



In this Section, we introduce the tight-binding model 
of the ABA-stacked trilayer graphene and briefly review 
the derivation of low-energy effective Hamiltonian in the 
absence of bias following Ref. M We will discuss that in 
the absence of perpendicular electric field (also referred to 
as displacement field), the band structure can be viewed 
as two independent bands: one bilaror-graphene-like and 
the other monolayer-graphene-like.'^SESnlSl 
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FIG. 1. (Color online) Schematic view of the lattice struc- 
ture of the AB A-stacked trilayer graphene (left) and values of 
corresponding tight-binding parameters adopted from Ref. [5] 
(right). Red (blue) atoms belong to A (B) sublattice of the 
corresponding layer. 



A. Tight-binding model 

In order to describe the band structure of the unbiased 
trilayer graphene, we adopt the standard Slonczewski- 
Weiss-McClure parametrization of the tight-binding 
model.'^ Corresponding tight-binding parameters are 
schematically illustrated in Fig. [T] Six tight-binding pa- 
rameters 7o • . . 75 describe hopping matrix elements be- 
tween different atoms: 



Ai^ Bi-. 70, 
^1 ^ ^3 : ^72, 

^1,3 ^ M . 
^1,3 ^ B2 



-81,30^2:71, (la) 
^1.3 O B2 : 73, (lb) 

Bi^B^: ^75, (Ic) 



where Ai (Bi) refers to an atom from A {B) sublattice, 
and index i = 1 . . . 3 labels three layers [see Fig. [1] for no- 
tation of sublattices] . In addition, parameter 5 accounts 
for an extra on-site potential for Bi, A2, and B^ sites, 
which are on top of each other. The values of differ- 
ent tight-binding parameters are far from being settled. 
In particular, authors of Ref. jSj used experimentally de- 
termined higher LLs crossings to refine values of tight 
binding parameters. Corresponding values are listed in 
Fig. [1] and will be used in what follows. 

Effect of an external perpendicular electric field (dis- 
placement field) can be described by including addi tional 
on-site potentials for different layers C/i . . . C/3 .mUEll, Jt 
is convenient to introduce parameters 



Ai = (-e 



A2 = (-e) 



Ui - 2U2 + 
6 



(2) 



Physically, parameter Ai is responsible for the potential 
difference between top and bottom layers which is caused 
by the application of the displacement field to the sam- 
ple. Parameter A2 describes the deviation of potential 
on the middle layer from the mean of the potentials on 
top and bottom layers. As we shall discuss below, non- 
zero parameter A2 can emerge when screening is taken 
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into account and the potential distribution is calculated 
self-consistently. In principle, even without bias, non-zero 
value of A2 is allowed by symmetry. In what follows we 
will demonstrate that LL crossings pattern is very sensi- 
tive to precise value of A2. 

For our analysis below, it will be convenient to separate 
the Hamiltonian of the unbiased trilayer, Hq, from terms 
related to the potential imbalance between layers. 



H — Hq + i/Ai + 



(3) 



In the momentum representation, the tight-binding 
Hamiltonian can be written in a compact matrix form, 
with different components corresponding to the wave 
function amplitudes on the six sublattices. Let us choose 
the basis {Ai, Bi, A2, B2, A3, ^3}. In this basis, the three 
terms in the Hamiltonian can be written as 



/ 7oifc litl l3tk 72/2 \ 

Jotk S 71 74tfc 75/2 

rr litk 71 S lQt*k li^k 71 

73^fc lAtk laik 73ifc ^^tk 

72/2 -iit*k l3tk 7ot; 

V 75/2 71 74t^ 7oifc S J 



Ha, =diag(Ai,Ai,0,0,-Ai,-Ai) 



(4) 



(5) 



i?A2 = diag(A2, A2, -2A2, -2A2, A2, A2). (6) 

In Eq. Q we used tk as a short-hand notation for a sum 
over nearest neighbors in a single honeycomb lattice: 



tk 



3 



(7) 



where ai = (0, l/\/3), 02,3 — (=Fl/2, — 1/2a/3), are vec- 
tors connecting an Ai site to its nearest neighbors, and 
momenta are measured in units of inverse lattice con- 
stant, a = 2.46 1. 

To obtain an effective low-energy Hamiltonian, we ex- 
pand Eq. Q in the vicinity of and if_ points, lo- 
cated at position (±47r/3,0) in the hexagonal Brillouin 
zone. This expansion reduces to replacing 7^^^. — UjTr, 
where 



TT = C^x + iky, hv., = 
with ^ = ±1 for and K_ points, respectively. 



(8) 



B. Effective Hamiltonian of unbiased trilayer 
graphene 



In the absence of external bias, the band structure 
of the Hamiltonian Q consists of monolayer- like and 
bilayer-like bands, as shown by Koshino and McCann.S^l 



To illustrate the decoupling of the Hamiltonian, we con- 
sider a different basis. 



^/2 



V2 



n/2 



n/2 



In the new basis, sum of two terms, Hq -\- H/\2 (A2 is 
allowed by symmetry even in the absence of electric field) , 
acquires a block-diagonal structure. 



Ho + Ha, 



( HsLG 

V 



H 



(10) 



BLG , 



where the monolayer-like and bilayer-like blocks are de- 
fined as 




The effective Hamiltonian of bilayer-like sector can be 
simplified further by noting that the low-energy states 
predominantly reside on Ai + A3 and on B2 sublat- 
tices. The higher energy bands, which reside mostly on 
sublattices Bi + B^, A2 have a large band gap of order 
±\/27i ~ 0.5 eV. Taking into account virtual excitations 
to these bands, one obtains an effective 2x2 Hamiltonian 
of the bilayer-like bandl^Sl 



H 



BLG 



r(0) 
'BLG 



^(1) 
'BLG' 



where h'^^q is the leading approximation. 



H 



(0) 
BLG 




1 

2m 



(12) 



(13) 



and ^^BLQ describes various corrections: 



-"blg ~ 



271 



/2i;3 



TT 

Trt 0, 



(5-2A2)7rt7r 




^ / 72/2 -f A2 
V -2A2^ 

° 1 

(75/2 + J + A2)wt^ 
r(i) 



(14) 



Note that in the above expression for -ffgLQ we omitted 

from terms which renormalize mass in h'^^q^ [Eq. (13)]. 
Also, the last term in Eq. (14) docs not include contri- 



butions which come with an additional small parameter 
7174/(70'^) « 0.11. 

Therefore, in the limit when there is no electric 
field (Ai = 0), band structure of ABA stacked trilayer 
graphene has SLG and BLG sectors which are decoupled 
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(a) Ai = (b) Ai = 0.025 eV (c) Aa = 0.05 eV (g) Ai = 0.2 eV 




FIG. 2. (Color online) Evolution of the band structure of trilayer graphene with bias Ai. The bands in the vicinity of 
A"+ point in the Brillouin zone are shown as a function of ka, where k describes momentum relative to the K+ point. Solid 
(dashed) lines correspond to k parallel to a;-axis (y-axis). At zero bias (panel (a)), red and blue lines describe monolayer-like 
and bilayer-like bands. At non-zero bias the monolayer-like and bilayer-like bands hybridize and no clear distinction can be 
made between them. As a result of this hybridization, one pair of bands moves to higher energies. The other pair of bands 
undergoes a series of transformations which lead to seven emergent Dirac points at higher bias Ai > 0.1 meV. The positions 
of these Dirac points in the Brillouin zone are illustrated in panel (g). 



from each other. The resulting band structure of unbi- 
ased trilayer is illustrated in Fig.[2]ja). 

When a displacement field is applied, the reflection 
symmetry between the top and bottom layer is broken 
and the H^^ appears in the Hamiltonian (note that elec- 
tric field also affects A2, which, in general, should be 
calculated self-consistently, see below). In the rotated 
basis ([9| this term corresponds to hybridization between 
the monolayer-like and bilayer-like blocks and has a form: 

/ i/hyb\ „ /Ai \ 

"""^^Wv. j' ^'''"lo aJ- 

(15) 

As a next step, we will discuss the evolution of the bands 
resulting from this hybridization. 



III. BAND STRUCTURE EVOLUTION UNDER 
BIAS: A NEW SET OF DIRAC POINTS 

An electric field completely changes the low-energy 
properties of trilayer graphene, hybridizing the mono- 
layer and bilayer-like bands. In this Section we shall 
study how the bands transform as bias increases. We 
will see that at relatively small bias, there is a complex 
interplay of terms induced by electric field and smaller 
trigonal warping terms in the Hamiltonian, which gives 
rise to rather flat low-energy bands. At higher bias 
(Ai > O.leV), a new set of Dirac points emerges: there 



is a total of seven DPs in each valley, one of them situ- 
ated at K± , and two groups of three related by rotational 
symmetry. The positions, mass gaps, as well as velocities 
of these DPs change as a function of bias. 

Since analytical treatment is no longer possible, we 
resort to numerical diagonalization of the tight-binding 
Hamiltonian ([4]). The evolution of the band structure 
in valley Kj^ as Ai is changed from zero to 0.25 eV is 
illustrated in Fig. [2j Notice that the band structure in 
K- valley is related to that in valley by the time- 
reversal symmetry. Throughout this Section, we will put 
A2 = 0, because A2 is expected to be of the order of 
several meV, and we found that it cannot qualitatively 
change the band structure. We emphasize, however, that 
A2 will play an important role for LL spectrum and LL 
crossing pattern in Section [IV[ 

Band structure in the absence of bias. Fig. [2ja), con- 
sists of a monolayer-like band (shown in red) and a trig- 
onally warped bilayer-like band (shown in blue). Rela- 
tively small bias Ai w 0.025 eV in Fig. ^h), turns out 
to be sufficient to hybridize bands; an electron-like and a 
hole-like bands which used to be of monolayer type, drift 
away to higher energies. The band gap for these bands 
becomes of the order 0.03 eV already at Ai = 0.025 eV, 
which means that they will be essentially irrelevant for 
the low-energy properties. 

As the bias is increased further (Ai — 0.05 eV in 
Fig.[2jc)), the two remaining bands stay close to the neu- 
trality point, and they are separated by a small band gap 
of the order of a few tens of meV. Near the iiT-points, the 
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bands have a complicated dispersion relation and trigo- 
nal warping is important. One notable feature is that 
the hole-like band is nearly flat in the vicinity of the 
point. 

When Ai becomes larger than 0.1 cV, the evolution 
of the bands leads to seven new Dirac points: one at the 
point, and six situated off-center. The off-center DPs 
consist of two groups of three, related by three-fold rota- 
tional symmetry. All Dirac points are generally split by a 
small mass. The locations and splittings of the emergent 
DPs for Ai = 0.2 eV are illustrated in Fig. ^g). In the 
interval of values Ai — 0.15-^ 0.25 eV, the off-center DPs 
move further away from the point. Also, for one 
of the groups the band gap closes and the Dirac mass 
changes sign at Ai w O.lSeV [see Fig. [2] (d)-(g)]. The 
velocity also changes, although not very significantly. 

To understand the origin of the new DPs and some 
of their properties, including their position in the Bril- 
louin zone, let us re-inspect the tight binding model at 
large bias. We note that at large bias, the most impor- 
tant terms in the band structure are 70, 7i, 73 and Ai. It 
is convenient to separate the Hamiltonian into the lead- 
ing part which contains only these terms, H^^\ and the 
remaining subleading part, H'^^^: 



(16) 



To elucidate the structure of ij'"^ let us re-arrange the 
basis ([9]) as follows: 



^1 -^^3 

V2 



B1+B3 A1+A3 



V2 



Bi — B3 



(17) 

Using Eqs. (|4|-([5|, we write in the vicinity of 
point in the new basis: 



C 

,ct 0, 



c = 



Ai 

VqTT 



VqTT' 

Voir 
V271 



(18) 

Similarly, one can find an explicit form of subleading part 
of H'-^\ but we will not need it here. 



The chiral symmetry of the Hamiltonian ( 18 ) greatly 



simplifies the analysis of its band structure. This symme- 
try implies that the spectrum is particle-hole-symmetric. 
We find that for any Ai, i?'") has seven massless DPs. 
Their positions can be found from condition det C = 0, 
which reduces to an algebraic equation on tt: 



AIvqtt + wq'"'* (27i'U37r* — VqTt'^) = 0. 



(19) 



First of all, tt = is an obvious root of the above equa- 
tion. Second, introducing parametrization tt = pe**, we 
obtain the other roots which give the positions of the 
remaining six off-center DPs: 



V7M + A|.§ +7i^^ ^ = 0,2vr/3,4,r/3,(20) 
"0 



P2 ^ ^ ' = 7r/3,7r, 57r/3. (21) 
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FIG. 3. (Color online) Evolution of the mass of emergent 
Dirac fermions as a function of bias for central DP (purple) , 
first group of off-center DPs (dotted green) and second group 
of ofi^-center DPs (dashed green). The mass of the second 
group of DPs vanishes at Ai « 0.185 eV and becomes negative 
at a larger bias. 



Therefore, at small Ai <^ JiVs/vq « 0.04eV the off- 
center DPs remain relatively close to the point 
and their position depends approximately quadratically 
on Ai: pi « 27iV3/wg + Al/{2-fiV3), P2 ~ A2/(27iW3). 
At large bias, Ai ^ 0.04 eV, the distance of the off- 
center DPs from grows linearly with Ai, Pi(2) ~ 
(Aiwo ±7i?;3)/wg. 

The chiral symmetry is broken by the subleading part 
of the Hamiltonian, ij'^^ in Eq. ( |l6[ ). Consequently, all 
seven species of emergent Dirac fermions become mas- 
sive, as is evident from Fig. [2] We extract the Dirac 
mass as a function of bias for these DPs from numerically 
evaluated spectrum. The result is illustrated in Fig. [s] 
The most interesting feature is that the mass of the sec- 
ond group of Dirac fermions vanishes at AI « 0.185eV, 
and changes sign for larger values of Ai. Note, that 
the value of A^ depends on the tight-binding parame- 
ters. The mass of the first group of Dirac fermions grows 
monotonically, changing from SmeV at Ai = 0.1 eV to 
ISmeV at Ai — 0.25 eV. Finally, the mass of the central 
DP depends on bias relatively weakly. 

The tunability of the Dirac masses will manifest it- 
self in a characteristic non-monotonic bias dependence of 
conductivity at the neutrality point. At Ai < A^ and 
Ai > A^ we expect an activated temperature depen- 
dence of the conductivity with an activation gap set by 
the smaller of the Dirac masses. When Ai « A^;, as- 
suming the disorder is weak, the conductivity will have a 
metallic temperature dependence and will be dominated 
by the twelve species of massless Dirac fermions (includ- 
ing spin and valley degeneracy), becoming of the order 
of 12 X e^/h (we assumed that a single Dirac specie has 
a conductivity approximately equal to e'^/fP). 

We note t hat, as previously pointed out by McCann 
and Koshino,'^l23 the screening can be significant, and 
therefore the parameters Ai,A2 should be calculated 
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self-consistently. We have calculated the self-consistent 
values of Ai, A2 as a function of external bias in the 
Hartree approximation, obtaining results which were con- 
sistent with Refs. |ni[27]. Interestingly, there is non-zero 
A2 ~ 2 meV even at zero external bias. This term will 
be important in our analysis of low-lying LLs in the fol- 
lowing Section. 



IV. LANDAU LEVELS AND THEIR 
CROSSINGS 

In this Section, we analyze Landau level spectrum of 
the biased trilayer graphene. As we shall see, the bias- 
induced band structure transformations lead to a rich 
pattern of LL crossings. At high bias, when the new 
Dirac points are formed, corresponding three-fold degen- 
erate groups of LLs emerge. This gives rise to an anoma- 
lous sequence of quantum Hall plateaus. At smaller bias, 
all LL degeneracies are lifted, and therefore all quantum 
Hall plateaus are spaced by /h. An interesting feature 
of the LL crossings pattern is its extreme sensitivity to 
certain tight-binding parameters, including A2 and 73. 
We explore how the existence and positions of various 
crossings depend on the values of these parameters, and 
argue that in an experiment the LLs crossing pattern can 
be used to determine the precise values of tight-binding 
parameters. 



A. Zero Landau levels without displacement field 

We start with reviewing known results about LL spec- 
trum in the absence of bias.--"' ------- In what follows, we 

shall be mostly interested in the low-lying LLs. The en- 
ergies of su ch LL s in the unbiased case can be obtained 
analyticalljiSHm from the effective low-energy model de- 
scribed in Section III AI This will facilitate the identifi- 
cation of various low-lying LLs in numerical simulations, 
and also will help us to develop an intuition about their 
sensitivity to tight-binding parameters. 



In the absence of bias, the monolayer-like [Eq. (llal] 



and bilayer-like [Eqs. (13l-(14)] bands are independent 
To find the spectrum of LLs, we adopt the Landau gauge. 
Ax = 0, Ay = Bx, and make a substitution tt — ^ tt — 
e{Ax + iAy) in Eqs. (11a) and (13)-(14). Momentum fc„ 



is conserved (and is proportional to the guiding center 
position X), and 7r-operator acquires the following form 
= ±1 for K± valleys): 



ih 



i^dx + x-x), 



(22) 



where the coordinate and momentum are measured in 
units of magnetic length Ib = ^Jh/eB and 1//b, respec- 
tively. We see that tt and tt^ can be viewed as lowering 
and raising operators of the magnetic oscillator. In the 
basis of non-relativistic LL orbitals \n) (at a given X), 



the matrix elements of tt, tt''' in the two valleys are given 
by, 



TT \n) 
TT^ \n) 
TT \n) 
TT^ \n) 



l-B 



+ l)|n + l) 



ih /— . 

— — V 2n |n 
Ib 

ih r-— , 

— V 2n \n — 



Ib 



(23a) 



(23b) 



The Hamiltonian of the monolayer-like band. 



Eq. (llal, describes two-dimensional Dirac- like fermions 



with a mass. Although the entire spectrum of Landau 
levels can be found analytically for the monolayer block, 
for our purposes it suffices to consider only the zeroth 
LLs in K± valleys: 



*0" 



£;o» =(5- ^75 + A2, (24a) 



(24b) 



Therefore, zeroth LL in monolayer sector is valley split, 
with the energy difference given by A™ = Eqi^ — Eq-^ = 
(5 — 75/2-I-72/2 w 7 meV (where we used tight-binding 
parameters listed in Fig. [T]). This situation should be 
contrasted with monolayer and bilayer graphene, where 
LLs are valley-degenerate, at least at the single-particle 
level. ^ 

Finding low-lying LLs in the bilayer block is more chal- 
lenging. In general, the trigonal warping described by the 
first term in ( 14 1 does not allow to solve for LL spectrum 



analytically, as it breaks rotational symmetry and ad- 
mixes arbitrarily high LLs. However, it turns out that 
trigonal warping gives relatively small corrections to LL 
energies, at least when magnetic field is not too small 
Therefore, one can find LL spectrum in the bilayer block 
neglecting trigonal warping, which gives 




, Eq^^ = -2A2, (25a) 
E,^^ = -2A2 + C (y + (J + A2):25b) 
+ A2, (25c) 
+ A2 + C('5-2A2), (25d) 



72 

2 

72 



where uJc — eB/m = ^/2eBv'^ /^i is the cyclotron fre- 
quency, and C — hujc/ {V^li) is a dimensionless param- 
eter proportional to the magnetic field [C, ~ 0.02 for 
B — 5 Tesla) that controls splitting between two LLs 
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(a) Ai = (b) S = 5 T (c) B = 10 T 




fi,T Ai,eV Ai,eV 



FIG. 4. (Color online) Landau level spectrum as a function of magnetic field and bias. Plot (a) shows LL spectrum as 
a function of magnetic field B — 2 ... 10 Tesla at zero displacement field, Ai — 0, and A2 = 0. Landau levels from the 
monolayer-like and bilayer-like block are shown in red and blue color; solid (dashed) lines correspond to K+ {K-) valleys. LLs 
are labelled with their indices, e.g. 1+, 1+ 1™ in the main text) stands for LL with index 1 from K+ valley of bilayer 
and monolayer sector respectively. Lifting of spin degeneracy by the Zeeman field is not shown, so that every LL is doubly 
degenerate. Neutrality point is located between 0+ and 1+ LLs. Plots (b)-(c) show Landau levels at _B = 5 T and i? = 10 T as 
a function of displacement field. Red (blue) color is used to indicate LLs which belong to monolayer (bilayer) sector at Ai — 0. 



in the same valley. The inter-vaUey sphtting, approxi- 
mately given by 72/2 + 3A2, is of the order of lOmeV; it 
is very sensitive to the value of A2, since it enters with a 
factor of 3. Below we will see that even relatively small 
variations of A2 can completely change the pattern of 
LL crossings at small negative filling factors. Also, no- 
tice that the inter-valley splitting is much larger than the 
intra-valley splitting: the latter is proportional to C, and 
is of the order 1 meV at B = 5 T. 



B. Landau levels at high bias 

We now proceed with the analysis of LL spectrum at 
a non-zero bias. As bias admixes the monolayer-like and 
bilayer-like blocks, and trigonal warping becomes impor- 
tant, approximate analytical treatment is no longer pos- 
sible. Therefore, we find LL spectrum and wave functions 
numerically. 

If we fix guiding center position at some value, our 
Hamiltonian can be viewed as an infinite matrix. It is eas- 
ily written in the basis of non-relativistic LL orbitals \n) , 
using matrix elements of tt, tt^ operators from Eq. (23). 



corresponding to raising and lowering operators accord- 
ing to Eq. (23 1 with the cutoff A. All other non-zero 



We truncate this infinite matrix, restricting the Hilbert 
space to LL orbitals with indices n < A. Finally, the 
problem of finding LL spectrum reduces to diagonaliz- 
ing a matrix in which each element of 6 x 6 Hamilto- 
nian, Eq. (10), is replaced by a matrix of dimensions 
A X A. Momentum operators are replaced by matrices 



elements 7^ are replaced by 7ilAxA- 

One subtlety associated with this procedure is due to 
the fact that the cutoff imposed on the ladder operators 
gives rise to new unphysical eigenvalues.'^fl' Such unphys- 
ical eigenvalues has to be identified and removed. We do 
this by using the fact that, although these "false" Lan- 
dau levels have low energy, their wave functions are dom- 
inated by Landau level orbitals with large indices (near 
cutoff A) . We have found that cutoff A > 30 is sufficient 
to faithfully represent the evolution of LLs in vicinity of 
neutrality point with magnetic field and bias. All simu- 
lation presented in this paper are done with A = 50. 

The evolution of LLs with magnetic field at zero bias 
is shown in Fig. |4] (a). The zeroth LLs in bilayer and 
monolayer block can be easily identified in this picture 
because they disperse very weakly with magnetic field. 
These results are consistent with Ref. [28], however, no- 
tice, that we have chosen a different set of tight-binding 
parameters extracted from experiment^, which leads to 
visible differences in the LL spectrum. 

Landau level spectrum as a function of bias for fixed 
values of magnetic field, B = 5T and i? = lOT, illus- 
trated in Fig. |4|^b)-(c), reveals a very complex pattern 
of crossings. At smaller field B — 5T and large bias 
A > 0.15 eV, both the dispersion and degeneracy of LLs 
changes: in particular, there are three-fold-degenerate 
groups which are associated with the new emergent Dirac 



points in the band structure. In addition, there is a 
singly-degenerate zeroth LL that originates from the cen- 
tral Dirac point in each valley. The dispersion of new 
zeroth LLs is determined by the dependence of the Dirac 
mass on the bias (see Fig. |3|. 

It should be noted that the three-fold degeneracy of 
a subset of DPs gives rise to the three-fold degener- 
ate groups of LLs only when magnetic field is not too 
strong. This means that the inverse magnetic length 
must be smaller than the distance between the DPs in 
the BZ. Considering Ai ^ 0.04 eV and using the ap- 
proximate expressions for the DP positions obtained in 
the previous Section, we write this condition as follows, 
h/ls ^ Ai/'^Q. Thus, three-fold degenerate groups of 
LLs are most easily observable at smaller fields. 

Experimentally, the degeneracy of certain LLs at high 
bias (and relatively small field) will manifest itself in an 
anomalous quantum Hall sequence. Some of the steps 
will be characterized by a jump of Hall conductivity by 
3e^/h (when a three-fold-degenerate LL is filled), while 
others will correspond to a jump of e^//i (when a singly 
degenerate LL is filled). Here we assumed that the spin 
degeneracy is lifted by the Zeeman interaction. In exper- 
iment, such an anomalous quantum Hall sequence will 
signal the formation of new DPs. 
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(a) 73 = 0.39 eV (b) 73 = 0.78 eV 
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FIG. 5. (Color online) Influence of the value of 73 on the 
pattern of avoided crossings between monolayer and bilayer- 
like LL. Magnetic field is B = 10 T. For clarity, only LLs in 
the K+ valley are shown. Plot (a) is for 73 = 0.39 eV used 
in this paper, whereas plot (b) is for twice larger value of 73. 
Labelling scheme of LLs is explained in the caption of Fig. |4] 



C. Anti-crossings and the role of 73 

As we discussed above, bias gives rise to three-fold de- 
generate groups of LLs, which is in full agreement with 
the emergencece of low-energy DPs in the band struc- 
ture. Another notable feature of a LLs pattern, which is 
visible in Fig. |4] (b)-(c) is that LLs from monolayer-like 
sector move away from neutrality point as the bias is in- 
creased. This fact can also be understood using intuition 
developed about the band structure evolution under the 
bias: hybridization and floating away of monolayer-like 
bands [see Fig. [2] (a)-(b)] corresponds to monolayer-like 
LLs moving away and exhibiting avoided crossings with 
bilayer-like LLs. 

Interestingly, it is tight binding parameter 73 that con- 
trols the behavior of aforementioned avoided crossings. 
For example, let us zoom in into part of LLs spectrum 
corresponding to negative filling factors of order v « —20. 
We note that such filling factors are easily accessible in 
the experiment.'^ Fig. [5] illustrates the pattern of avoided 
crossings for moderate bias Ai = 0.17 eV and two 
different values of 73: the one which is used in all simu- 
lations throughout this paper [73 = 0.39 eV, Fig. [s] (a)] 
and a twice larger value [73 = 0.78 eV, Fig. [s] (b)]. 

In Fig. [sjjb), corresponding to the larger value of 73, 
the crossing pattern changes due to a slightly different 
order of LLs at zero bias. More noteworthy is the pro- 
nounced enhancement of gaps at anti-crossings. Due 
to this enhancement, one can identify the clear pattern 
in Fig. [5] (b): avoided crossings occur between a given 
LL originating from monolayer-like band and every third 



bilayer-like LL. This pattern can be easily explained us- 
ing the effect of trigonal warping on the LL wave func- 
tions. If, without a trigonal warping term, a given LL 
contains orbital \i) in its wave function, a non-zero trig- 
onal warping admixes orbitals \i ± 3), |i ± 6), . . . with 
progressively smaller coefficients. It is this admixture 
that underlies the avoided crossings, thus, it is natural 
that stronger trigonal warping enhances the gaps at anti- 
crossings. 

As an example, let us consider LL 1™ in the monolayer- 
like block [labeling scheme of LLs is introduced in the 
caption of Fig. |4] . From Fig. [5] we see that it exhibits 
a pronounced anti-crossing with 8*^^ bilayer-like LL, and 
also with ll^J. LL for larger 73. In the absence of trigonal 
warping, the wave function of 1™ LL at Ai = consists 
of a combination of |0) and |1) orbitals on Ai — A3 and 
Bi — B3 sublattices respectively [see Eq. (|9|]. The wave 
function of 8*^ LL is dominated by |6) and 78) orbitals on 
symmetric combination of A1.3 sublattices and B2 sublat- 
tice. Non-zero trigonal warping changes the wave func- 
tions, admixing substantial components of |0), |3), |9), . . . 
orbitals to the |6) orbital on Ai + A3 sublattice. Fi- 
nally, non-zero Ai breaks inversion symmetry, and ini- 
tially orthogonal symmetric and antisymmetric combi- 
nation of Ai,3 sublattices evolve to have non-vanishing 
overlap. This gives rise to the anti-crossings observed in 
Fig. [5] 

The observed sensitivity of the pattern of avoided 
crossings can be potentially used to fix the value of the 
tight-binding parameter 73. However, in order to do this. 
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FIG. 6. (Color online) Landau levels as a function of displacement field for different values of A2. Magnetic field i? = 10 T. 
Plot (a) uses zero value of A2, for plots (b) and (c) A2 = 3 meV, and A2 — 7 meV respectively. Value of A2 = 7 meV is 
sufficient to switch the order of (1,0)_ and (1,0)+ LLs from bilayer sector. 



the other tight-binding parameters have to be refined. 
Below we discuss how information about them can be 
extracted from the LLs crossings at small filling factors 
and a smaller bias. 



D. Pattern of Landau level crossings at smaller 
bias and its sensitivity to tight-binding parameters 

Having discussed the behavior of LL at high bias and 
(or) large filling factors, we concentrate on the regime of 
small bias, when the new Dirac points have not formed 
yet. To better understand the complicated pattern of 
LL crossings that arises in this regime, we zoom in on 
the region Ai < O.leV, and concentrate only on several 
low-lying LLs [see Fig. [6] (a)]. 

Interestingly, the bias completely lifts the degeneracies 
of all LLs. Therefore, we expect that in biased bilayer 
graphene all quantum Hall plateaus will be spaced by 
e'^/h. Notice that this goes against intuition acquired 
from the studies of monolayer and bilayer graphene,!^! 
where in the absence of interactions Landau levels re- 
main valley-degenerate. 

Landau levels undergo a series of crossings, as illus- 
trated in Fig. [6|a). In particular, level 0*1 moves down, 
crossing weakly split 2+ LLs at Ai w 0.03 eV. The en- 
ergy of level 1*1 increases with bias, such that it crosses 
level 1*^ at Ai ss 0.07 eV. Furthermore, levels 0^ and 1^ 
cross at Ai « 0.05 eV as 1+ level slowly moves down, and 
0+ level is essentially non-dispersive. Finally, monolayer- 
like levels 0™ and 0+ both float up and cross at small 



bias Ai w 0.015 eV. They also cross 2+ LLs which move 
down. 

As we emphasized above, although qualitatively the 
band structure of trilayer graphene is rather well un- 
derstood at this point, not all parameters in the tight- 
binding model are precisely known. The complex LL 
crossing pattern described above is very sensitive to the 
choice of tight-binding parameters, and may provide a 
valuable tool for fixing their values. We illustrate this 
by exploring how LL crossing pattern changes when pa- 
rameter A2 is varied. Such a parameter is allowed by 
symmetry, but have not been taken into account previ- 
ously. 

The LL crossing pattern for A2 = 3meV and A2 = 
7meV is illustrated in Fig. [6]jb)-(c). It is evident that 
the crossing between l^i and 1+ is significantly shifted 
even by small A2 = 3 meV; so does the crossing between 
monolayer-like 0+ LL and 2+ levels. 

As A2 is increased further [see Fig.[6](c)], new crossings 
appear: in particular, a pair 0+, 1+ swaps positions with 
the pair 0*1, 1*1 at zero bias [as expected from Eq. (25)]. 



This leads to new crossings between 0_ and 0+, 1+ levels 
at a small bias, Ai w 0.01 eV. Also, 0+ LL moves above 
2± levels, such that crossings between them disappear. 
We see that even small variations of parameter A2 can 
change the pattern of LL crossings. Similarly, the posi- 
tions of different crossings are sensitive to the values of 
72, 75, S, although in a more subtle manner. 

As we pointed out, non-zero A2 can be present among 
other tight-binding parameters. Moreove, if one accounts 
for the effect of electric field self-consistently, in addition 
to renormalization of Ai (corresponding to screening). 
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parameter A2 changes as well. We found that typical 
values of A2 obtained from such self-consistent procedure 
is 2meV at charge neutrality point We also performed 
numerical simulations for the case when magnetic field is 
present (detailed procedure will be described elsewhere). 
We found that the self-consistent value of A2 in the pres- 
ence of magnetic field can be enhanced compared to the 
case when B = 0, due to discreteness of LLs. For ex- 



ample, from Eqs. (25a)- (25b I one can see that at zero 



bias wave functions of bilayer-like 0!j. and 1^ LLs are lo- 
calized on the middle layer. Consequently, by emptying 
these LLs one can create an excess positive charge on the 
middle layer, thus generating positive non-zero A2 that 
can be as large as 4meV.'^ 



V. DISCUSSION AND OUTLOOK 

In summary, we have studied the electronic properties 
of the AB A-sta.cked trilayer graphene in the perpendicu- 
lar electric field. We found that the hybridization of the 
monolayer-like and bilayer-like bands leads to an opening 
of a small (< 10 meV) band gap at experimentally achiev- 
able values of bias Ai < 0.25 eV. At bias Ai > 0.1 eV, 
a new set of Dirac points appears in each valley, which 
includes one central DP, and six off-center DPs. Masses 
of the emergent Dirac fermions are tunable by bias, and, 
interestingly, some of the masses vanish at a bias which 
for our choice of tight-binding parameters corresponds 
to AJ « 0.185 eV. Therefore, the band gap depends on 
the bias in a non-monotonic manner and closes at the 
bias A|. This behavior should be contrasted with bi- 
layer graphene,^ where band gap grows monotonically 
as a function of bias, and can become very large, of the 
order of 0.25 eV.El 

We have also studied the spectrum of LLs, finding that 
at bias Ai > 0.1 eV the band structure transformation 
gives rise to new three-fold-degenerate groups of LLs. 
The evolution of low-lying LLs from the unbiased case, 
where the spectrum consists of LLs in the monolayer and 
bilayer sectors, to the degenerate groups of LLs at higher 
bias, is accompanied by multiple level crossings. Experi- 
mentally, these LL crossings will give rise to phase transi- 
tions between quantum Hall states as a function of bias. 
We have also shown that the pattern of crossings depends 
strongly on the values of certain tight-binding parame- 
ters, such as A2 and 73. Studying the phase transitions 
in QHE regime experimentally should enable the precise 
determination of these parameters. 

Our results show that single-particle splittings between 
LLs in biased trilayer are sizable, of the order of a few 
meV, and exceed the Zeeman energy of 1 meV at typi- 
cal magnetic field of i? = 10 T. This is expected to play 
an important role in the analysis of qu antum Hall fer- 
romagnetism in trilayer graphene,'22Ji2il which will be a 
subject of subsequent publication.'^ Here we just note 



that we found that, owing to the small LL spacing, the 
effective Coulomb interactions are strongly screened,!^ 
similar to the case of bilayer graphene^^ and, as a re- 
sult, the Coulomb energy scale becomes comparable to 
typical single-particle energy splitting. This situation 
is quite different from the case of monolayer graphene, 
where single-particle splittings of LLs are much smaller 
than the interaction energy. ■^^ •^^ 

Our study reveals the unique tunability of Landau 
levels in trilayer graphene: bias controls their energies, 
wave functions, and degeneracies. We speculate that 
this may enable the exploration of interesting interaction- 
induced phenomena in the QHE regime. First, the three- 
fold degeneracy of LLs at higher bias may allow one 
to study quantum Hall ferromagnetism, skyrmions, and 
fractional quantum Hall states with S'[/(3)-symmetry. 
Second, the bands in trilayer graphene are strongly trig- 
onally warped, which provides an opportunity to study 
the effect of anisotropic band structure on quantum Hall 
ferromagnetism and fractional QHE.I^SI Third, Landau 
level wave functions involve a superposition of different 
non-relativistic LL orbitals with weights tunable by bias. 
Therefore the LL form-factors are tunable, which may al- 
low one to realize new effective interaction regimes, sta- 
bilize desired fractional quantu m Ha ll states and drive 
phase transitions between them.'2212l] 

We believe that phenomena discussed above can be 
experimentally observed in the near future. We note 
that the observation of effects associated with the emer- 
gent DPs is possible when disorder broadening is smaller 
than the energy scale over which effective description in 
terms of DPs is valid. Fig. [2] indicates that the latter 
energy scale is of the order of 20meV, which is much 
larger than the disorder broadening achievable in trilayer 
graphene on hexagonal boron nitride.'^'^ Therefore, new 
DPs should be observable in currently available samples. 

Finally, we note that very recently multiple phase tran- 
sitions between different quantum Hall states in biased 
trilayer graphene have been reported.!^ It is possible that 
these phase transitions can be understood in terms of 
single-particle LL spectrum evolution described in our 
paper. 
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